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Abstract 

The KdV equation can be considered as a special case of the general equation 

ut + f{u)^-5g{uxx)x = 0, 6 > 0, (0.1) 

where / is non-linear and g is linear, namely f{u) = and g{v) = v. As the 
parameter 6 tends to 0, the dispersive behavior of the KdV equation has been 



throughly investigated (see, e.g., ||ll|, [y], |^ and the references therein). We 
show through numerical evidence that a completely different, dissipative behavior 
occurs when g is non-linear, namely when g is an even concave function such as 
g{v) = —\v\ or g{v) = —v"^. In particular, our numerical results hint that as (5 — > 
the solutions strongly converge to the unique entropy solution of the formal limit 
equation, in total contrast with the solutions of the KdV equation. 
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1 Introduction 

Let us consider the first-order non-linear evolution PDE 

+ = 0, (1.1) 

with a strictly convex non-linearity f"{u) > a > 0, typically f{u) = v? /2 which cor- 
responds to the usual inviscid Burgers equation |TT|. This equation generally produces 



shock waves in a finite time and solutions must be understood in a suitable weak sense. 
A good framework enforcing existence, uniqueness, stability, (|l^, [§, ID,.--) is 
provided by the so-called entropy solutions, obtained through the vanishing viscosity 
method by passing to the limit in the viscous approximation 

ut + /(m)x - = 0, (1.2) 
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as (5 > tends to zero |T0[. Entropy solutions can be also characterized as the weak 
solutions of which satisfy the Oleinik one-sided Lipschitz condition (OSLC) 

Ux{t,x) < (1.3) 
at 

A completely different behavior is observed if ( |1.1| ) is approximated by the KdV equation 
ut + f{u)x - Suxxx = 0, 5 >0, (1.4) 

as investigated by many authors, in particular ||^. When the formal limit has shocks, 
the corresponding solution of the KdV equation becomes highly oscillatory and does not 
converge at all to the entropy solution (even in a weak sense). 

An alternative approach of approximating equations such as ( p..lD is to design a 
numerical scheme. In [|1| Brenier and Osher presented a second-order accurate version 
of the Roe scheme {^^, |§]) and an OSLC consistent method of lines for approximating 
solutions of scalar conservation laws, ( |1.1|) . As was already pointed out in the 
"modified equation" of this scheme (that is to say the asymptotic PDE of the discrete 
equation) is of interest since (as it is second-order) it must involve, just like the KdV 
equation, a third-order term, but unlike the KdV equation, due to the enforced OSLC 
this third-order term must enforce restrictions on the possible creation of wiggles. Our 
present work presents an extended numerical study of these "non-classical dispersive" 
effects, as well as an extended study of the numerics involved. 

The model equations we consider are ( p.l| ) where g is concave and even, typically 
g{s) = —\s\ and the smoother g[s) = — s^. The numerical tools taken from ||l[| enable 
us to approximate the solutions of ( p.l| ) with a numerical method in which its modified 
equation takes the desired form of 

Ut + f{u)x - 5g{uxx)x = -euxxxx + LOT, (1.5) 

with a positive e ~ 6 Ax, where LOT means perturbation terms with at most third-order 
space derivatives (low order terms). In the specific case where g{s) = — \s\, there are no 
low order terms on the RHS of ( p..5|) and hence we are left only with a linear dissipative 
term proportional to 6 Ax, which enables us to hold a fixed 6 while eliminating the RHS 
at the limit e — > 0. 

The analytical tools available to study fully non-liner evolution PDEs such as ( p.l[ ) 
and ( |1.5| ), seem to be limited and a rigorous theory is, in our opinion, a very challenging 
task that deserves to be further addressed. Our main results are based on numerical 
evidence. 

The simple non-linear switch in the third-order term lu^xU based on the sign of 
the second derivative Uxx, has a major impact. The same holds true for other non- 
linear third-order operators such as {ul^)x- Instead of the usual dispersive wave-train 
that develops in the KdV equation, here we still observe the appearance of a single 
ripple, which indeed can be still considered as a dispersive phenomenon. This ripple, 
however, is of a totally different character than what we are familiar with in the KdV 
equation. Instead of developing into a full wave-train which propagates in time, this 
ripple remains single during the evolution of the solution. It does not increase above a 
certain size which is related to the value of the parameter 6 in (|0.1|). 
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The striking numerical discovery is that not only this ripple is single, but it also 
disappears in the limit 5 — > 0. It turns that (p.l|) is much closer in its behavior to the 
dissipative, viscous approximations of the form 

ut + f{u)x = eu^x, e > 0, 
than to the KdV equation 

Ut + f{u)x + Su^xx = 0, 

for which at the same limit 6^0, the dispersive waves increase. 

The natural conjecture that follows is that ( |1.5|) , strongly converges to the entropy 
solution of the scalar conservation law (|1 . 1|) as 6 and e6~^ tend to 0. 

We would like to emphasize that when using numerical methods in order to approx- 
imate solutions of differential equations, instead of solving the original problem, we end 
up solving a modified problem, which typically includes high-order terms. Here, we 
are interested in studying the effects of such high-order terms in the particular case in 
hand. We believe that the properties of the model equations we study in this work, 
might provide insight for the inclusion of similar nonlinear terms in future models. 

The paper is organized as follows: we start in §^ by presenting both our model 
equations and the numerical method of approximating their solutions. We continue 
in ^ by commenting on the major properties of these equations. In particular, we 
study traveling waves solutions and state a formal one-sided bound (OSLC-type) on the 
derivative. 

We then continue in §^ by presenting several numerical examples which leads us to 
our main conjecture stated above, which is formalized in Conjecture [4.1| . 

We end in §^ by demonstrating and analyzing spurious numerical oscillations that 
appear when solving our model problems with a certain choice of parameters. The 
methods we use for analyzing the binary oscillations follow Goodman and Lax, [§. 
Similar techniques were used, e.g., in ^, |^. This study is aimed at gaining a further 
understanding in order to be able to distinguish between the properties of the equations 
and the appearing numerical phenomena which depend on the specific method used 
when approximating the solutions. 

Acknowledgment: The research was supported in part by TMR grant ^^ERBFM- 
RXCT960033. The work of D.L. was supported in part by the Applied Mathematical 
Sciences subprogram of the Office of Science, U.S. Department of Energy, under con- 
tract DE-AC03-76-SF00098. Part of the research was done while D.L. was affiliated 
with ENS Paris. D.L. would like to thank E. Tadmor and A. Kurganov for stimulating 
discussions. 

2 Motivation 

2.1 The Model Equations 

The motivation for our study is the modified equation of the Brenier-Osher (BO) scheme 
presented in We start by considering the scalar conservation law (|1 . 1|) with a 
(strictly) convex flux function, f"{u) > a > 0. 
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We discretize the space using a fixed mesh spacing, Ax, denoting the approximation 
of the solution at the discrete grid-points by := u{iAx). 

The BO method for approximating solutions of ( |1.1D is a second-order convergent 
Roe scheme and an OSLO consistent method of lines, which can be written as 

^ + ^(/.«/.-/-V.)=0, (2.1) 

with 

/i+l/2 — n \Ui -\ —SiyUi+i 2""^*^^ 

Here, h is the Godunov flux associated with /, 

h{a, h) = haooio,, b) = e min e/(a -|- s{b — a)), e = sgn{b — a), 

0<s<l 

and the slopes, s,, are limited by 

.. = ^max(«.^,-«„«.-«._0. 

We note that the OSLO consistency of (|2.1| ) is in the sense that the discrete approx- 
imation satisfies 

Ui+i{t) - Ui{t) ^ 1 



Ax - p(m(0, + 
where p{v) is the one-sided Lipschitz semi-norm, 

/ vix) — viy) \ 
p[v) = sup I , 0+ = max(0,a). 



x^y \ X y 

In the linear case, i.e. for f{u) = (3u with (3 constant, the BO scheme yields the 
modified equation 

ut + f{u)x + 5\uxx\x = Q, 5 = 0(Ax2). (2.2) 



In the nonlinear case for a general /(m), additional low-order terms appear in ( |2.2| ). The 
existence of these terms was ignored in |Q. 

As was already pointed out in [Q, this modified equation ( |2.2| ) is of interest since it 
involves a dispersion term (as it is second-order), but unlike KdV, due to the enforced 
OSLO, this third-order term does not create wiggles. These are the properties this work 
is aimed at studying. 

Since ( p.2|) includes a non-smooth flux, we extend the current discussion to a wider 
family 

ut + f{u)x - 5g{uxx)x = 0, (2.3) 

where giuxx) is assumed to be concave, g"{s) < 0,Vs. There are no further regularity 
assumptions on g. 

To summarize, our model equations are ( p.2| ) and its smooth version ( p.3| ). 
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2.2 The Numerical Scheme 



For the following analysis the precise derivation and the origin of ( |2.2| ) and ( p.3| ) is no 
longer important. It is their properties that we want to study. To proceed, we write 



a discretization of ( p.2| ),( |2y3| ). It is remarkable that it is possible to write down such a 
discretization whose modified equation involved a linear fourth-order derivative of the 
"correct" (dissipative) sign. 

Following IQ, a straightforward discretization of the generalized modified equation 
( ^731) will be, e.g., to define the (first-order) numerical flux gi+1/2 = g{wi,Wi+i) as the 
Engquist-Osher (EO) flux 



gEo{a,b) = g{b^) + g{a_), 6+ = max(0,6), a_ = min(0,a). 
A method of lines for (E.3|) can be then written as 



(2.4) 



at 



f 



6g{wi,w. 



0, 



(2.5) 



where A_aj = ai — aj_i, and Wi = (ui_i — 2ui + Uj+i)/(Ax)^. 

Replacing the time derivative in ( p. 51 ) by the forward-Euler discretization associated 
with a fixed time-step. At, we obtain the fully-discrete 



H A 

At Ax 



f 



(2.6) 



where for the space-time discretization we use the notation m" := u{iAx, nAt). 

A somewhat tedious but straightforward computation shows that the limiting equa- 
tion for the fully-discrete process ( |2.6|) has the form of 



ut + f{u)x - 5g{Uxx)x 
= Ax6 



(2.7) 



^sgn{uxx)g' {uxx)uxxxx + ^sgn{uxx)g"{uxx)u. 



2 

XXX 



+ 0{{Axf,At). 



In particular, for the choice g{s) = — (|2.7| ) becomes 

1 



™, e-^Ax5. 



Ut + f{u)x + 5\Uxx\x = -eUx 

For the more regular choice of g{s) = — s^, (|2.7D takes the form of 



(2.8) 



Ut + f{u)x + 5{ul^)x = -Ax5 



\Uxx\Uxxxx "t~ ^^Sg^iUxxjU^ 



(2.9) 



At that point, we add equations (|2.8| ),( pl9D to our two previous model equations, 
(23) and dg^). 
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3 Properties of the Model Equation 

In this section we discuss some of the major properties of (|2.3| ). We will assume that 
f{u) is convex. For demonstration purposes if will be usually taken as f{u) = u'^/2. 
We also assume a concave g{s), which will be either — |s| or — s^, unless mentioned 
otherwise. Namely, we deal with two model equations. The first equation corresponds 
to the choice of g{s) = — \s\, 

Ut + UU^ + 6\Uxx\:^ = 0. (3.1) 

The second choice corresponds to the smoother, g{s) = — s^, 

ut + uux + 6{ul^)x = 0. (3.2) 

We note that in both equations ( |3.1D and (|3.2|), one can eliminate by rescaling the 
parameter 6. This is done in ( p.l|) by the rescaling u —>■ 6u and t —>■ 6t. Such a rescaling 
also holds for ( p^.8| ) when e ~ 5 Ax. The different scaling x S^^^^x and t 5^^^^t gets 
rid of 5 in (|33). 

We would also like to note that Ut + uUx + Sg{uxx)x = —^"^xxxx is invariant under the 
Galilean transformation ^(a;, t) — >■ u{x + ct, t) and u ^ u — c. 



3.1 Traveling Waves 

We compute a monotonic profile of the traveling waves solution for ( ^.2| ). The mono- 
tonicity is an important property of these traveling waves. Our numerical experiments 
below, indicate that these particular traveling wave solutions behave as attractors. 

We impose a symmetric upstream-downstream profile with Ui^^a = — Mmght = ui > 0. 
One integration of (|3.2|) results with 



- - ^ = -6u' 
2 2 

which can be integrated twice to yield 

-1/2 



X 



■(2(5)^/4 £ (^s^ul -s^ + ul arcsin (^^^ - "^u]^ ^ ds. (3.3) 



The sign in the integrand is set according to the sign of u. In particular, for mi = 1, 
rSD becomes 



X 



■(2(5)^/^^ (^sv^r^ + arcsin s- I)) ^'"^ ds. (3.4) 



We demonstrate this profile in the following example in which we approximate so- 
lutions of (|3.2|) using the numerical scheme presented in § |2.2| . Here, we used Riemann 
initial data as 

%(^) = \ _\ (3.5) 
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The results are plotted in figures |3Tl|^l2| for different values of 6, after the steady- 
state was reached. These numerical results hint that the traveling wave solutions behave 
as attractors to the approximate solutions. 



numerical solution 
traveling wave 



Figure 3.1: Convergence of the Numerical Solution to the Traveling Wave Solution, 5 

10-5, T = 0.5 



1 'K^^^-o--<^^^-o--e-^^^^^^^^e^«^>^o-,o. 



numerical solution 
traveling wave 



-0.3 -0.2 



0.1 0.2 



Figure 3.2: Convergence of the Numerical Solution to the Traveling Wave Solution, S 

10-6, T = 0.5 
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3.2 The One-Sided Lipschitz Condition 

It is surprisingly easy to derive a formal one-sided bound on the derivative Ux in (|2.3| ). 
Here, we have to augment ( |2.3| ) with additional assumptions on g{s) given in the fol- 
lowing elementary lemma 

Lemma 3.1 (OSLC) Consider ^7^ ). Assume thatWs, f"{s) > 0,g"{s) < and that 
g'{0) = 0. Then the derivative, u^, satisfies the formal one-sided estimate 

max-u^(t) < maxM2.(0), Vt > 0. (3.6) 

X X 

Proof: Differentiating ( p.3| ) once yields 

utx + {f'{u)ux)x - S{g'{uxx)uxxx)x = 0. (3.7) 

Denoting w = Ux, { \i-7\) can be rewritten as 

wt + f"iu)w^ + f'{u)wx - 5g"{wx)wl^ - 5g'{wx)wxxx = 0. 

If we now take v{t) = max^,. w(a;, t), then due to the assumptions of the lemma, v < 
and ( |3.6| ) follows. ■ 

The appearance of a fourth-order derivative in ( |2.7| )-( PTP| ) prevents us from writing 
an analog of lemma in this case. The fact that the fourth-order derivative appears in 
the "correct" negative sign in the RHS of the equations is associated with a dissipative 
process, but to quantify this in terms of the desired bound on the derivative seems to 
be a highly non-trivial task. 



4 Strong Convergence? 



In this section we numerically study the limit of solutions of ( p.2|) -(2.3) as the parameter 
S tends to zero. Our starting point is the results obtained by approximating solutions of 
( p.l|) subject to the exponentially decaying initial conditions u{x, 0) = e-i°°^ . Similar 
phenomena to those we describe below, appear with other choices of initial data. 

Figure [4.1| presents our most interesting discovery. One can observe the creation 
of a single dispersive "ripple" which disappears as 5 — > 0. This phenomenon is totally 
different from what we know about solutions of dispersive equations (KdV-type). There, 
an infinite dispersive wave-train develops in time (compare, e.g., with p. 64]). This 
wave-train increases at the same limit 5 — > 0. In that limit, the solution appears to 
strongly converge to the entropy solution of the Burgers equation subject to the same 
initial data. 



We confirm that the results shown in figure |4.1| are independent of the spatial spacing 
Ax, by demonstrating in figure O a case in which 6 is being held fixed while approaching 
the limit Ax 0. 



The evolution in time of the solution is shown in figure We note that the 

dispersive ripple does not increase in magnitude above a certain size (which by figure 
is determined by the value of 6). No other ripples are being developed. 
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In figure we present a fully numerical phenomenon in which spurious oscillations 
propagate from the shock. This happens when the equation is numerically solved for a 
relatively small value of S and a relatively large Ax. We observe that these oscillations 
disappear as As 0. A further study of similar numerical phenomena is presented in 

m 

A single ripple which disappears at 6 tends to zero is observed when equation ( |3.1j ) is 
replaced by the smoother variant (|3.2| ) subject to the same exponential initial data. An 
example of the approximated solutions of ( |3.2| ) is presented in figure Once again we 
observe the appearance of a single dispersive ripple which vanishes in the limit 6^0. 
As before, in that limit, the solution appears to strongly converge to the unique entropy 
solution of the Burgers equation subject to the same initial data. 

Our numerical results in this section lead us to the following 



Conjecture 4.1 The solution of 



ut + f{u)x - Sg{ua:x)x = -euxxxx, 5, e > 0, 



with a strictly convex f and a concave g, with g{0) = 0, strongly converges to the unique 
entropy solution of as 6 and e6~^ tend to 0. 





exact 
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0.001 
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0.6 



0.5 
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0.2 



0.1 



-0.1 



-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 0.6 O.i 



Figure 4.1: Absolute- Value Dispersion with Exponential Initial Data. = 400, T = 0.5, 
At ■ 5 = 10~^. The different plots refer to different values of 5. The "exact" curve is the 
entropy solution of the Burgers equation subject to the same initial data 



10 



Y. Brenier and D. Levy 





Figure 4.3: Absolute-Value Dispersion with Exponential Initial Data. Time evolution of the 
solution. N = 400, 5 = 5- 10"^, At ^2- 10"^ 
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Figure 4.4: Absolute- Value Dispersion with Exponential Initial Data. The appearance of 
numerical oscillations for small 6 = 10"^ T = 0.5. For N = 200, At = 10"^. For 
N = 100, At = 10-3 
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0.0001 
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Figure 4.5: Smooth Dispersion with Exponential Initial Data. = 400, T = 0.5, At ■ 5 = 
10~*. The different plots refer to different values of S The "exact" curve is the entropy 
solution of the Burgers equation subject to the same initial data 



5 Study of the Spurious Numerical Oscillations 

In this section we demonstrate and analyze several cases in which numerical oscillations 
appear. We emphasize that the goal of this study is to gain a better understanding in 
order to distinguish between the properties of the model equations and the phenomena 
which are strictly related to the numerics in hand. 

In figures pTTHsT^ we present results obtained by approximating solutions of ( p. 3D 
with g{uxx) = —\uxx\, subject to the Riemann initial data ( |3.5| ) for several different sets 
of parameters. 

One can clearly observe the binary oscillations propagating from the discontinuity at 
the center of the domain. These oscillations are strictly associated with the numerical 
method. 

Note, in particular, the difference between the two oscillatory solutions presented in 
figure ^.2| . While the support of the oscillations of the solution for = 100 expand in 
time, the solution for N = 200 is a steady-state. Even though we do not know what 
mechanism toggles between these two states, we are able in the steady-state case to 
compute the envelope of the decaying oscillations. 

We would like to stress that these oscillations appear in the non-smooth g{s) = —\s\, 
but do not appear with the smoother g{s) = — s^. Such oscillations-free results were 
already presented in the study of traveling waves with the smooth g{s) = — s^, (Figures 
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The numerical initial data which was used in figures |5ri| - |5-3| includes the origin (0, 0). 
The existence of this point in the initial data has a strong stabilizing eff'ect on the results. 
It seems to control the oscillations such that they do not change sign and they remain 
bounded in time. This is not the case when this extra point is not used. For comparison 
we bring in figures |5.4| - |5.(j| equivalent results to figures [5.1| - p.3| when no such extra point 
is used. Similar sensitivity of the oscillations with respect to the location of the grid 
points was previously observed in a different context by Hayes in p|. 




-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 0.6 0, 



Figure 5.1: Absolute- Value Dispersion with Riemann Initial Data. S = 10 ^, T = 0.5 
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Figure 5.3: Absolute-Value Dispersion with Riemann Initial Data. S — 10 ^, N — 100 
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-1 -0.8 -0.6 -0.4 -0.2 



Figure 5.4: Absolute- Value Dispersion with Riemann Initial Data which does not include the 
origin. 5 = 10-^ T = 0.5 




Figure 5.5: Absolute-Value Dispersion with Riemann Initial Data which does not include the 
origin. 5 = 10-^ T = 0.5 
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Figure 5.6: Absolute- Value Dispersion with Riemann Initial Data which does not include the 
origin. 5 = 10-^ N = 100 



We turn to a study of the binary oscillations presented in the above examples. The 
method we use follows the lines of Goodman and Lax in |^. Similar techniques were 
used, e.g., in §1, 0]. 

Our starting point is the method of lines (|2.5|) with the specific choice of f{u) = u'^/2. 
We assume that the numerical solution introduces binary oscillations of the form 



Ui = Vi + {-lywi. 

We plug ( p.l| ) into ( |2.5| ). For the first term, we have 



Ax 



-A 



/ 



8Ax 



[{Ui+l + Ui-l + 2Mj)(Mi+i - 



Vj+i - Wj+i + Vj^i - Wj-i + 2vi + 2wi Vj+i - Vj-i - Wj+i + Wj-i 
4 2Ax 

Vj+i + Wj+i + Vj-i + Wj-i + 2vi - 2wi Vj+i - Vj^i + Wj+i - Wj-i 
4 2Ax 



(5.1) 



(5.2) 



I even 



i odd. 



If we consider equation ( ^.2|) as a discrete approximation of continuous equation, it 
can be rewritten as 



Ax 



-A_ 



/ 



Ui + u 



i+l 



v[Vx — Wx) t even 
v{vx + Wx) i odd. 



(5.3) 
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For the second term in ( p.5| ) we repeat the same analysis and end with 

A^g{wi,Wi+i) — > 



(5.4) 



9 V, 



Aw 



(Ax) 



2 ) '^xx 



9 [v^x + 



Aw 
(Ax) 



Aw 
(A^ 

Aw 
(Ax)2 



5- ^^zx + 



- 5- - 



4w 



(Ax) 



Aw 

Aw 
(Ax)2 



« even 



i odd. 



When we specifically choose the fiux-limiter as the Engquist-Osher (EO) limiter, 
becomes 



A_g{wi,Wi+i) 



sgn{a)g{a) — sgn{[3)g{[3) i even 
—sgn{a)g{a) + sgn{[3)g{[3) i odd, 



where 



a = Vxx + 



Aw 
{Axj 



/5 = f XX - 



Aw 
(A^ 



(5.5) 



(5.6) 



We require that the continuous solution will satisfy both equation for the even and odd 
indices. Combining the results of ( [5. 3D and ( |5.5| ), the continuous equations that corre- 
spond to the discretization of the method of lines ( p.5|) given by the binary oscillations 
solution (15. 11) are 



{V + W)t + V{Vr, 



Ax 



{sgn{a)g{a) - sgn{j3)g{j3)) = 0, 



(5.7) 



{v - w)t + v{v^ + w^) - —{-sgn{a)g{a) + sgn{f3)g{f3)) = 0, 



where a and (3 are given in ( ^.6|) . If we specifically chose g{s) 
4w/(Ax)^, ( |5.7| ) becomes 



\s\ then for Vxx ^ 



[v + w)t + v{vx - Wx) + 

{v - w)t + v{Vx + Wx) - 



8(5 



(Ax)3 
(A^) 



w = Q, 



(5.^ 



= 0. 



Here, it is clear that by taking the limit 5 — > 0, in order not to affect the oscillations 
that might develop, the ratio (5/(Ax)^ has to be held fixed. Of course, numerically in 
this case, the RHS has a term —euxxxx-, with e = Ax5 which is of order (5/(Ax)^. Taking 
5 large enough for a fixed Ax means that there is more dissipation and in such a way 
oscillations can be controlled. 

An example of such sort of control of the oscillations is shown in figure |5.7| . Here, 
we approximated the solution of (|3.1|) , subject to the Riemann initial data ( ^.51) . The 
ratio 5/(Ax)^ was held fixed. 
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Figure 5.7: Steady-State Solution, T — 0.5. Fixed 5/{Ax)^ Controls the Oscillations 



As previously noted, in several cases when binary oscillations appear, it seems as 
if there is convergence to steady-state solutions. The envelope of this steady-state 
solutions can be computed from the above analysis, since what we are left with in 
this case amounts to (assuming positive w > 0) 



Vt + vv^ = 0, 

Wt — VWx + cw — 0, 

where c = 8S/{Ax)^. Then, for dt — 0, we have 

v'v'^ = 0, 

-v^wt. + cw^ = 0. 



(5.9) 



(5.10) 



Hence, = ±c which can be determined by taking the limit x ±oo. For the 
particular choice oi — ±1 one has 

We experimentally observed that when the initial data includes the origin, the oscilla- 
tions do not change sign and hence c is limited by the size of the discontinuity. The 
envelope of the oscillations is therefore given by 



V ±w = ±l±e (-^^'^ 



Fully Non-Linear KdV-Type Equations 



19 



where the different variations of the sign correspond to the upper /lower envelope at each 
part of the solution. 

In figures [578|-|0| we plot the results obtained by solving (|0| ) subject to the Riemann 



initial data ( |3.5| ). The numerical solution convergence to an oscillatory steady-state 
solution in which the oscillation are binary. 

In figure p.8| , the parameters used were 6 = 0.00001, = 100. The corresponding 
outer envelope takes the form 



-1 + e"''^ X > 
1 — e""^^ X < 0, 

while the inner envelope is 

—1 — e~^^ X > 
1 + e'"^ X <0. 



c= 10, 



c= 10, 



In figure |5.9| , the parameters were changed to 6 = 0.000001, = 200. The corre- 
sponding outer envelope takes the same form, with a different constant c = 8. 



1 1 1 1 J 


steady state numerical solution 




inner envelope 






outer envelope 













-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 0.6 0.8 1 



Figure 5.8: The Envelope of the Steady-State Oscillatory Solution, N = 100, 6 = 10 ^, 
T = 0.5 
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Figure 5.9: The Envelope of the Steady-State Oscillatory Solution, N — 200, 5 — 10 ^, 
T = 0.5 
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